# ==========================================
# Inter-rater Reliability Analysis in R
# Fleiss' Kappa and Pairwise Cohen's Kappa
# ==========================================

# Install required packages
install.packages(c("readxl", "irr", "dplyr", "openxlsx"))

# Load packages
library(readxl)
library(irr)
library(dplyr)
library(openxlsx)

# 1. Import dataset
data <- read_excel("Fleiss_Kappa_Jamovi.xlsx")

# 2. Check column names
names(data)

# 3. Select reviewer classifications
ratings <- data %>%
  select(A, B, C, D)

# 4. Convert ratings to factors
ratings <- as.data.frame(lapply(ratings, as.factor))

# 5. Overall Fleiss' Kappa
fleiss <- kappam.fleiss(ratings)

fleiss_results <- data.frame(
  Method = "Fleiss Kappa",
  N = nrow(ratings),
  Reviewers = ncol(ratings),
  Kappa = fleiss$value,
  Z = fleiss$statistic,
  p_value = fleiss$p.value
)

print(fleiss_results)

# 6. Pairwise Cohen's Kappa between reviewers
reviewer_pairs <- combn(names(ratings), 2, simplify = FALSE)

pairwise_kappa <- lapply(reviewer_pairs, function(pair) {

  result <- kappa2(ratings[, pair])

  data.frame(
    Reviewers = paste(pair, collapse = " vs "),
    Kappa = result$value,
    Z = result$statistic,
    p_value = result$p.value
  )

}) %>%
  bind_rows()

print(pairwise_kappa)

# 7. Simple percentage agreement
full_agreement <- apply(
  ratings,
  1,
  function(x) length(unique(x)) == 1
)

agreement_results <- data.frame(
  N = nrow(ratings),
  Reviewers = ncol(ratings),
  Agreement_Percentage = mean(full_agreement) * 100,
  Cases_Full_Agreement = sum(full_agreement),
  Cases_Disagreement = sum(!full_agreement)
)

print(agreement_results)

# 8. Identify records with disagreement
disagreement_cases <- data %>%
  mutate(
    Full_Agreement = full_agreement
  ) %>%
  filter(Full_Agreement == FALSE)

# 9. Pairwise contingency tables
pairwise_tables <- lapply(reviewer_pairs, function(pair) {

  tbl <- table(
    ratings[[pair[1]]],
    ratings[[pair[2]]]
  )

  as.data.frame.matrix(tbl)

})

names(pairwise_tables) <- sapply(
  reviewer_pairs,
  function(pair) paste(pair, collapse = "_vs_")
)

# 10. Export results to Excel
wb <- createWorkbook()

addWorksheet(wb, "Fleiss_Kappa")
writeData(wb, "Fleiss_Kappa", fleiss_results)

addWorksheet(wb, "Pairwise_Kappa")
writeData(wb, "Pairwise_Kappa", pairwise_kappa)

addWorksheet(wb, "Agreement")
writeData(wb, "Agreement", agreement_results)

addWorksheet(wb, "Disagreement_Cases")
writeData(wb, "Disagreement_Cases", disagreement_cases)

for (sheet_name in names(pairwise_tables)) {

  addWorksheet(wb, sheet_name)

  writeData(
    wb,
    sheet_name,
    pairwise_tables[[sheet_name]],
    rowNames = TRUE
  )
}

saveWorkbook(
  wb,
  "Interrater_Reliability_Results.xlsx",
  overwrite = TRUE
)

# 11. Summary output
cat("\nFile exported: Interrater_Reliability_Results.xlsx\n")
cat("Fleiss' Kappa:", round(fleiss$value, 3), "\n")
cat("Simple agreement:", round(agreement_results$Agreement_Percentage, 1), "%\n")
cat("Disagreement cases:", agreement_results$Cases_Disagreement, "\n")